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SECTION  I 


SUMMARY 

This  report  reviews  the  following  four  tasks  pertaining  to 
response  equalization  problems: 

•  Evaluating  a  new  technique  using  large  signals 
for  equalizing  seismometers 

•  Examining  the  concept  of  designing  regional 
equalization  filters  for  subarray  outputs 

•  Analyzing  statistically  the  coefficients  used 
to  equalize  the  noise  data 

•  Developing  the  theory  to  incorporate  statistical 
phase  fluctuations  in  the  correlation  statistics 

Section  II  discusses  the  technique  using  large  signals  to 
equalize  seismometers.  The  method  is  based  on  finding  the  minimum-phase 
filters  which  equalize  the  signal  power  spectra  for  a  set  of  seismometers. 
The  method  is  tested  by  using  the  nine  in-line  seismometers  (NW-SE  arm) 
of  LASA  subarray  C2  for  a  large  Aleutian  Islands  event.  Four  sets  of  filters 
are  designed  (using  two  gate  lengths)  to  equalize  eight  channels  to  the  refer¬ 
ence  channel  (seismometer  10).  The  reference  channel  is  subtracted  from 
the  other  channels  before  and  after  equalization  and  the  error  traces  com¬ 
pared.  The  equalized  signals  are  significantly  more  similar  for  the  first 
few  cycles  of  the  arrival.  After  the  first  few  cycles,  there  is  little  differ¬ 
ence  between  the  two  sets  of  signals  —  probably  because  of  interfering 
scattered  energy.  The  technique  appears  to  be  valid,  but  a  complete  evalu¬ 
ation  (including  an  analysis  of  cost  vs  improvement)  is  recommended. 

Section  III  discusses  the  designing  of  a  set  of  regional  equali¬ 
zation  filters  for  LASA  subarray  outputs.  Average  Levinson  filters  are  de¬ 
signed  for  three  subarray  outputs,  using  three  events  from  the  Aleutian 
Islands  region.  Then,  the  filters  are  applied  to  seven  events  from  this 
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region,  and  results  are  compared  with  those  obtained  for  both  amplitude 
equalization  and  individual  Levinson  equalization.  Regional  equalization  fil¬ 
tering  is  found  to  be  possible  — -  but  only  if  the  epicentral  region  is  very 
small.  Thus,  regional  equalization  filtering  does  not  appear  to  be  practical 
to  implement. 

Section  IV  shows  that  variations  in  the  coefficients  used  to 
equalize  the  noise  are  caused  by  a  combination  of  seismometer  gain  fluctua¬ 
tions  and  variations  in  estimating  the  noise  average  power  (i.  e.  ,  the  zero- 
lag  autocorrelation  function  estimate).  Thus,  the  anticipated  use  of  the  noise 
equalization  coefficients  to  study  the  nature  of  statistical  gain  fluctuations 
is  not  possible. 

Section  V  shows  how  statistical  phase  fluctuations  can  be 
included  in  che  correlation  statistics.  This  is  a  generalization  of  gain  fluctua¬ 
tion  (  a  1 -point  stochastic  filter)  to  an  n-point  stochastic  filter.  If  the  filter 
weights  are  independent,  it  is  shown  that  phase  fluctuation  is  accounted  for 
in  the  same  manner  as  gain  fluctuation  (i.  e.  ,  scaling  the  autocorrelation 
function).  The  scalar  is  1  plus  the  sum  of  the  variances  of  the  n  points. 
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SECTION  II 


SEISMOMETER  EQUALIZATION  USING 
LARGE  SIGNALS 


A.  INTRODUCTION 

The  problem  of  equalizing  seismometers  occurs  frequently  in 
data  analysis;  for  example,  seismometer  equalization  is  necessary  in  comput¬ 
ing  wavenumber  spectral  estimates  and  in  designing  multichannel  filter  sys¬ 
tems  from  measured  noise  dat£  and  a  tneoretical  signal  model.  Techniques 
for  equalizing  seismometers  using  large  signals  range  from  simple  amplitude 
equalization  to  Levinson  equalization.  The  method  chosen  by  the  analyst 
depends  on  the  nature  of  his  problem;  each  method  offers  advantages  and  is 
subject  to  limitation^ . 

This  section  presents  a  new  method  using  large  signals  to 
equalize  seismometers.  Given  a  reference  channel  and  a  channel  to  be  equal¬ 
ized,  this  method  finds  the  minimum-phase  filter  which,  when  applied  to  the 
channel  to  be  equalized,  causes  it  to  have  the  same  power  spectrum  as  the 
reference  channel;  i.  e.  ,  if  gr  (t)  and  g  (t)  are  the  two  channels  with  power 

spectra  G  (f)  and  G  (f),  respectively,  the  desired  filter  is  the  minimum- 
r  r  6 

phase  filter  with  power  response  Gr  (f)/Ge  (*)• 

Differences  in  seismometer  responses  are  due  to  differences 
in  instrument  and  amplifier  responses  and  in  subsurface  structure.  Instrument 
response  and  impulse  response  of  a  layered  medium  are  minimum  phase. 
(Amplifier  response  may  not  be.)  Additionally,  the  difference  between  two 
minimum-phase  systems  is  minimum  phase,  so  the  choice  of  a  minimum- 
phase  filter  appears  to  be  a  reasonable  approximation.  Note  that  the  filter 
which  is  calculated  will  equalize  the  two  signals  only  if  the  difference  between 
them  is  actually  minimum  phase;  also,  the  two  signals  will  not  be  equalized 
if  the  spectral  ratio  changes  outside  the  design  gate. 
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A  program  to  design  and  apply  these  equalization  filters  uses 
a  new  technique  for  estimating  power  spectra.  Developed  by  Texas  Instru¬ 
ments  Incorporated,  it  reduces  end  effects  and  is  particularly  suited  to  short- 
duration  signals.  Filter  design  and  application  are  actually  accomplished 
entirely  in  the  time  domain  by  applying  to  the  channel  to  be  equalized  its 
own  whitening  filter  and  the  inverse  of  the  refe  re  nee -trace  whitening  filter. 

Figure  II- 1  shows  the  results  of  a  test  case  used  to  check 
the  method;  the  third  trace  is  the  actual  signal  arrival,  and  the  first  trace  is 
the  output  of  a  predetermined  minimum- phase  filter  applied  to  the  signal 
arrival.  Using  the  first  trace  as  the  reference,  the  third  trace  is  then 
equalized  using  the  equalization  program.  Note  that  a  very  short  gate  is 
used  *o  estimate  the  power  spectra  of  the  two  traces.  The  second  trace 
shows  the  estimate  obtained.  Considering  that  only  a  very  short  gate  and  a 
10-point  filter  are  used,  a  very  good  estimate  of  the  target  trace  was  obtained. 

B.  DATA  PRESENTATION 

Figure  U-2  shows  the  Andreanof  Islands  event  which  was  used 
to  evaluate  the  method.  The  nine  in-line  seismometers  of  subarray  C2  (Fig¬ 
ure  U-3)  were  chosen,  with  seismometer  10  as  the  reference  trace. 

Signal  sample  rate  is  0.  05  sec;  the  finer  sampling  is  preferred 
to  the  normal  0.  1  sec  to  improve  overall  accuracy,  especially  when  aligning 
the  signals.  Figure  II-4  shows  that  the  signal  peaks  slightly  above  1.  0  cps. 
The  signal  is  amplitude -equalized  (using  a  900-point  gate)  and  static-correc¬ 
ted  to  the  nearest  integer  sample  using  the  signal- reference  crosscorrelation 
function  (Figure  II- 5).  The  bottom  nine  traces  in  Figure  II-5  are  the  error 
traces  obtained  by  subtracting  the  reference  channel  (seismometer  10)  from 
the  other  channels.  (The  error  traces  are  scaled  up  by  2,  relative  to  the 
signal  traces.)  Note  that  the  error  traces  become  larger  as  seismometer 
separation  increases.  To  evaluate  the  equalization  method,  the  difference 
traces  are  compared  with  those  obtained  after  equalization. 
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Figure  II  - 1.  Results  of  a  Test  CaseUsed 
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.  Original  Data,  Andreanof  Islands  Event,  Subarray  C2 
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Figure  II-5.  Amplitude -Equalized,  Static-Corrected  Andreanof 
Islands  Event 
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To  make  valid  the  comparison  between  error-trace  sets 
it  is  necessary  to  align  the  signal  and  reference  traces  as  exactly  as  possible 
before  subtraction.  Assume  that  we  have  two  unit- amplitude,  monochromatic, 
1-cps  cosine  waves,  and  that  one  wave  is  displaced  by  A  sec  (A  <  At,  the  sam¬ 
ple  rate)  relative  o  the  other.  To  examine  the  effect  of  the  misalignment  on 
the  error  trace,  we  compare  the  average  power  in  the  error  trace  with  that 
in  the  original  trace.  Let 

Xj  (t)  =  cos  2  tt  t 

X2  (t)  =  cos  2  tt  (t+  A)  since  A  =  f  =  1 

Let  Y  (t)  be  the  difference  trace: 

Y  (t)  sX1(t)-X2(t) 

=  cos  2  tt  t  -  cos  2  tt  (t+  A) 
a  (COS  2  TT  t)  (  1  -  COS  2  TT  A) 

+  (sin  2TT  t)  (sin  2  TT  A) 

For  the  original  trace, 

X(t)2  =  f  cos2  2  TT  t  dt  since  T  =  1 

*0 

t_  sin  4  tt  t  I  * 

=  2  8TT  \z 

_  l 
~  2 
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For  the  error  trace 


Y  (t) 


=  f  icos  2 rr t ( 1  -  cos  2  tt  A) 

*'o  ' 

TT  A^ 


+  sin  2  TT t  sin  2 
1 


dt 


=  (1  -  cos  2  rr  A)2/  cos“  2  TT  tdt 
0 

2  Z*1  .  2  , 

+  sin  2  TTA  /  sm  ..  FT  t  dt 


i 

+  2  (1  -  cos  2  TTA)  sin  2TTA  f  cos  2tt  t  sin  2  TT  tdt 

.2  i  t  .  sin  4tt  t  I1  ) 

=  (1  -  cos  2  tt  A)  j-  +  8tt  I  j 


2  _  .  t  .  sin  4TTt 

+  sin  2  ttA'  j  +  — 3— 


=  i  |l+cos2  2  TTA  +  sin2  2 TTA  -  2  cos  2tta| 


=  -|  [2  -  2  cos  2TT  A] 


Thus, 

Y(t  f  /  X(t,2  =  2  ••  2  cos  2tt  A 

Assuming  that  A  is  one-half  block  (or,  in  our  case,  0.025  sec), 
the  error-trace  average  power  would  be  just  16  db  less  than  that  of  the  original 
trace.  The  observed  error  traces  are  6  db  to  18  db  below  the  signal  traces, 
so  a  one-half  block  misalignment  error  would  contribute  significantly  to  the 

error-trace  power. 
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Since  trace  alignment  is  critical,  a  program  has  been  written 
in  order  to  interpolate  the  signal  between  samples  so  that  the  maximum  mis- 
alignment  will  be  much  less  than  one-half  sample.  Appendix  A  outlines  the 
program  which  is  referred  to  as  the  microstatics  program.  Figure  11-6  shows 
that  the  error  traces  obta.ned  for  the  amplitude-equalised  data  are  smaller  - 
sometimes  by  a  significant  amount  (e.g.,  seismometers  46,  26,  and  53)  - 
after  application  of  microstatics. 

The  amplitude -equalized  microstatic-corrected  traces  of 
Figure  1-5  are  used  as  the  "original"  data  for  comparison  with  the  outputs 
of  the  equalization  filters  designed.  To  insure  a  valid  comparison,  the 
equalized  data  are  also  amplitude -equalized  and  microstaticly  corrected 
before  the  error  traces  are  obtained. 

Figures  11-7  through  H-9  show  the  equalized  traces  for  10-,  20- 
and  30-point  filters,  respectively.  Each  figure  shows  the  design  gate;  note 
that  the  same  signal  segment  is  used  in  each  case.  There  appears  to  be  no 
consistent  difference  in  the  performances  of  the  three  filters;  the  filter  length 
that  gives  the  smallest  error  trace  varies  from  channel  to  channel.  A  com¬ 
parison  of  the  equalized  and  original  error  traces  shows  that  equalization 
has  significantly  improved  the  signal  similarity  over  the  first  few  cycles; 
beyond  the  first  few  cycles,  there  is  little  difference  between  the  sets  of 
traces.  A  possible  explanation  for  this  behavior  is  that  scattered  energy 
accounts  for  a  large  proportion  of  the  difference  between  the  signal  traces; 
i.e.  ,  the  first  few  cycles  of  the  signal  arrival  are  free  of  scattered-energy, 
so  the  differences  are  largely  due  to  the  seismometer  responses,  which  the 
filters  correct  quite  well.  However,  the  differences  due  to  scattered  energy 
later  become  significant  with  respect  to  seismometer  response  differences; 
and  no  overall  improvement  is  obtained.  The  design  gate  used  to  determine 
the  filters  includes  much  more  of  the  signal  than  just  the  first  few  cycles. 
This  implies  that  the  filters  are  not  especially  attuned  to  the  first  few  cycles 
but  indeed  use  general  spectral  information  and  the  minimum-phase  assump¬ 
tion  to  achieve  equalization. 
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Figure  II-6.  Comparison  of  Error  Traces  for  Amplitude- 
Equalized  Data  Before  and  After  Microstatics 
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Figure  II-7.  Equalled  Signal  Using  a  10-Point  Filter 
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Figure  II-9.  Equalized  Signal  Using  a  30-Point  Filter 
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Figure  11-10  shows  the  equalized  and  error  traces  for  a  f>0- 
point  filter  designed  from  a  600-point  gate.  Again  equalization  is  most  ef¬ 
fective  for  the  first  few  cycles.  Also,  the  error  traces  are  slightly  lower 
overall  than  are  the  original  error  traces;  the  reason  for  this  is  not  known. 

Figure  II- 1 1  shows  the  signal-to-reference-trace  spectral 
ratios  before  and  after  equalization.  As  can  be  seen,  the  filters  equalize 
the  signal  spectra  quite  well.  The  large  positive  ratios  at  high  frequency 
before  equalization  are  due  to  the  lower  noise  *evel  for  the  deeper  reference 
trace  (seismometer  10).  Thus,  results  might  have  been  better  with  a  different 
seismometer  as  the  reference. 

C.  CONCLUSIONS 

Equalization  significantly  improves  signal  similarity  for 
the  first  few  cycles  of  the  signal  arrival.  Later  in  the  arrival,  no  improve¬ 
ment  is  obtained,  possibly  due  to  the  effect  of  scattered  energy. 

Since  the  equalization  technique  works  only  if  the  assumption 
that  signal  differences  are  minimum  phase  is  met,  it  appears  that  actual 
differences  are  approximately  minimum  phase.  Further  study,  therefore, 
is  recommended  to  more  fully  evaluate  the  method. 
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SECTION  III 


EVALUATION  OF  REGIONAL  EQUALIZATION 
FILTERS  vOR  LASA  SUBARRAY  OUTPUTS 

The  purpose  of  this  study  is  to  determine  the  feasibility  of 
designing  for  LASA  subarray  outputs  a  set  of  equalization  filters  which  can 
be  applied  to  all  events  from  a  given  epicentral  region.  If  regional  equaliza¬ 
tion  were  feasible,  a  set  of  equalization  filters  ior  each  region  could  be  de¬ 
signed  and  the  correct  set  automatically  applied  to  every  event  before  large- 
array  processing. 

A.  METHOD  OF  ANALYSIS 

Table  III- 1  shows  seven  events  from  the  Aleutian  Islands 
area  which  were  chosen  for  the  analysis.  Subarray  outputs  from  Bl,  B2, 

C3,  and  D2  (Figures  HI-1  through  IH-7)  were  used  for  each  event.  These 
subarrays  were  chosen  because  the  movecut  anomalies  were  the  same  for 
each  event,  thus  allowing  the  averaging  of  crosscorrelation  data.  The  events 
have  been  resampled  to  a  0.  1-sec  rate  and  bandpass -filtered  with  a  zero- 
phase  0.  8-  to  2.  8-cps  digital  filter  (Figure  IH- 8).  Table  IH-2  lists  the 
peak-signal-to-rms-noise  ratios. 

In  all  cases,  events  are  equalized  to  the  subarray  Bl  output. 
The  three  equalization  procedures  are  as  follows. 

•  Amplitude  Equalization 

The  square  root  of  the  ratio  of  the  energy  on 
the  reference  trace  to  that  on  the  trace  to  be 
equalized  is  computed,  and  the  trace  to  be 
<  dualized  is  scaled  by  the  resulting  number. 

•  Levinson  Equalization 

For  each  event,  21-point  Levinson  equalization 
filters  are  designed  using  100-point  (10- sec) 
gates.  The  filters  are  then  applied  to  the  ap¬ 
propriate  traces. 
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November  1965  Andreanof  Islands  Event 
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Figure  III- 8. 


Amplitude  Response  of  3.7-sec 
Filter  (100-msec  Sample  Rate; 


Zero-Phase  Bandpass 
0.  8-  to  2.  8-cps  Passband) 


•  Average  Levinson  Equalization 

Correlation  data  for  the  three  largest  events 
(17,  24,  and  36)  are  scaled  to  be  approximately 
the  same  size;  auto-  and  crosscorrelations 
are  stacked  and  an  "average"  Levinson  filter 
designed  for  each  of  the  B2,  C3,  and  D2  outputs. 

The  average  filters  are  then  applied  to  the 
traces  for  all  events. 

To  compare  Levinson  equalisation  with  amplitude  and  individual 
Levinson  equalization,  the  correlation  coefficient*  between  B1  and  the  other 
three  subarray  outputs  is  computed  for  each  set  of  equalized  traces.  Note 
that  scaling  does  not  affect  the  correlation  coefficient,  so  the  correlation  co¬ 
efficients  after  amplitude  equalization  are  the  same  as  the  correlation  coef- 
ficients  for  the  original  data. 
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B.  DATA  PRESENTATION 


Figured  Hi-9  through  III-  15  show  the  amplitude -equalized, 
the  Levinson-equalized,  and  the  average  Levinson-equalized  traces  for  each 
event.  Table  1H-3  lists  the  correlation  coefficients  obtained  after  each  type 
of  equalization.  The  events  are  grouped  by  epicentral  location;  the  first  four 
are  from  the  Andreanof  Islands,  and  the  last  three  are  from  other  locations 

in  the  Aleutians. 

A  comparison  of  the  raw  data  (Figures  HI-1  through  IH-7) 
with  the  amplitude -equalized  data  (top  set  of  traces  in  Figures  IH-9  through 
HI- 15)  shows  that  amplitude  equalization  is  necessary  for  all  events.  Ampli¬ 
tude  equalization  would  probably  be  sufficient  also  for  most  processing  schemes 
for  events  13,  17,  35,  and  36,  since  the  correlation  coefficients  are  close 

to  unity. 

Levinson  equalization  improves  the  signal  similarity  for  all 
events.  All  correlation  coefficients  except  those  for  subarray  D2,  event  24, 
and  subarray  03,  event  40,  are  greater  than  0.9  after  Levinson  equalization. 

The  increased  waveform  similarity  is  most  evident  for  event  40  (Figure  HI- 15). 

Average  Levinson  equalization  improves  the  signal  similarity 
for  events  13,  17,  24,  and  40,  does  not  appreciably  change  the  similarity  for 
event  35,  and  decreases  the  similarity  for  events  33  and  36.  For  the  four 
events  where  similarity  is  improved,  the  improvement  is  less  than  that  ob¬ 
tained  by  Levinson  equalization.  Note  that  these  four  events  have  almost 
identical  epicenters  (within  approximately  100-km  distance);  the  three  re¬ 
maining  events  have  significantly  different  epicenters  (up  to  700-km  distance). 

It  thus  appears  that  average  equalization  is  possible  but  that  average  filters 
are  critically  dependent  on  event  location  and  therefore  can  be  applied  only 
to  events  from  a  very  small  epicentral  region.  This  strong  dependence  of 
the  filter  set  on  event  location  appears  to  make  average  equalization  impractical 
because  of  the  large  number  of  filter  sets  that  would  be  required. 
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Figure  IU-9.  Equalized  Data,  Event  13 
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Table  III- 3 

CORRELATION  COEFFICIENTS  FOR 
ORIGINAL  AND  EQUALIZED  DATA 


Event 

No. 

Subarray 

A 

(°) 

Azimuth 

n 

Correlation  Coefficients 

Orig  &  Amp  Eq 

Lev  Eq 

Avg  Lev  Eq 

13 

B2 

46.  1 

304.0 

0.907 

0.975 

0.964 

C3 

0.852 

0.971 

0.  914 

D2 

0.930 

0.930 

0.  924 

17* 

32 

46.  5 

304.  1 

0.887 

0.975 

0.958 

C3 

0.  920 

0.950 

0.928 

D2 

0.  902 

0.945 

0.931 

& 

24 

B2 

45.  9 

303.6 

0.754 

0.931 

0.915 

C3 

0.722 

0.942 

0.884 

D2 

0.723 

0.803 

0.  772 

40 

B2 

46.8 

304.  5 

0.756 

0.969 

0.901 

C3 

0.763 

0.873 

0.846 

D2 

0.811 

0.948 

0.840 

33 

B2 

41.  2 

300.8 

0.851 

0.944 

0.803 

C3 

0.791 

0.936 

0.  ?5 1 

D2 

0.825 

0.917 

0.824 

35 

B2 

39.  2 

302.6 

0.  930 

0.991 

0.  928 

C3 

0.960 

0.992 

0.975 

D2 

0.  975 

0.991 

0.  975 

36* 

B2 

51.  4 

309.4 

0.943 

0.980 

0.871 

C3 

0.911 

0.949 

0.876 

D2 

0.849 

0.959 

0.836 

* 


Used  in  design  of  average  Levinson  filter 
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Note  that  ev?nt  36  has  poorer  signal  similarity  after  average 
equalization,  even  though  it  is  one  of  the  three  (with  events  17  and  24)  used 
in  the  average  filter  design.  Apparently,  events  17  and  24  together  domi¬ 
nate  in  the  average  filter  design.  It  appears  that,  if  event  36  had  been  omit¬ 
ted  from  the  average  filter  design,  correlation  coefficients  for  events  13,  17, 
24,  and  40  after  average  equalization  would  have  been  almost  as  good  as 
those  obtained  after  individual  equalization  for  these  events. 

C.  CONCLUSIONS 

* 

As  observed  previously,  amplitude  equalization  of  LASA 
subarray  outputs  would  generally  be  sufficient  for  most  processing  schemes. 

The  design  of  a  set  of  regional  equalization  filters  appears  to 
be  feasible  but  only  if  the  epicentral  region  is  very  small.  Consequently, 
implementation  seems  impractical  except  for  specific  small  areas  of  interest. 


SECTION  IV 

ANALYSIS  OF  NOISE  PREEQUALIZATION  COEFFICIENTS 


A  previous  report  shows  that  the  noise  at  a  subarray  must  be 

equalized  prior  to  multichannel  processing  in  order  to  obtain  consistent 

* 

rejection  at  low  frequency  (below  1.  0  cps).  The  technique  adjusts  the  rms 
noise  level  on  each  trace  to  a  common  value;  i.  e.  ,  the  quantity 


where 


a.  is  rms  noise  level  on  i  trace 

l 

.th  ,  ,th 

X..  is  l  sample  on  i  trace 
Ji 

N  is  number  of  samples  in  each  trace 


NC  is  number  of  channels  in  a  subarray 


is  computed  for  every  trace  and  the  quantity  is  found  such  that 


Pi 


a 

r 


J. 

l 


where 


a  and  a. 
r  l 


is  scalar  to  apply  to  i  trace 

are  rms  noise  levels  on  reference  trace  and  i 
trace,  respectively 
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This  technique  is  simple  to  implement  and,  because  of  the  peaked  noise 
power  spectra,  effectively  equalizes  the  data  at  low  frequency.  Consequently, 
it  was  applied  at  all  subarrays  for  all  noise  samples. 

The  large  number  of  equalization  coefficients  collected  appear 
to  provide  the  data  necessary  to  undertake  a  statistical  analysis  of  seismom¬ 
eter-gain  fluctuations.  Suppose  we  express  the  noise  on  any  channel  at  any 
time  as 


X..  =  K.  X.. 
Ji  i  Ji 


where 


K.  is  gain  of  ith  seismometer-amplifier  system 

.th 


X.'  is  ground  motion  at  i  seismometer  at  time  j 


K.  is  assumed  to  be  a  random  variable  which  is  statistically  independent  of 
the  gain  on  any  other  seismometer.  From  physical  considerations,  it  appears 
reasonable  to  assume  that  K.  is  constant  over  the  time  interval  (0,N);  i.  e.  , 
the  seismometer  gain  is  assumed  to  be  a  slowly  varying  function.  Then 


we  can  express  as 


Pi  ~  0 


N 


a  i  E  (k  X, 

r  J  N  1=0  v  r  jr/ 


1/2 


N 


(k.x./) 

N  .  n  '  i  Ji 7 
J=° 


K 

_ L  J 

irl 

i 

i  e  (xji')2  J 

1/2 


Kr  I  *rr  (°> 


1/2 


Ki  ]  fa  (0) 
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where 


are  zero-lag  values  of  autocorrelation 
function  on  channels  r  and  i,  respectively 


cp  (0)  and  cp. .  (0) 
rr  11 


or 


K 


P.  = 


i  K. 


a 


where 


a  =  [cPrr  (0)  j  ^(0)] 


1/2 


if  a  is  nearly  constant  for  each  noise  sample,  for  an  assumed  distrubution  of 
the  K's,  the  distribution  of  the  p's  can  be  calculated.  (For  example,  if  the 
K's  are  uniformly  distributed  anda«  1,  distribution  of  the  p's  can  be  derived 
as  in  Appendix  B. )  It  would  then  be  possible  to  compare  the  theoretical  and 
observed  distribution  of  the  p's. 

The  key  to  the  problem  is  to  determine  whether  the  assumption 
that  a  is  nearly  constant  is  reasonable;  estimating  the  variance  of  a  is  necessary. 
Because  the  time  series  is  autocorrelated  (i.  e.  ,  the  X..'s  are  not  independent), 
the  exact  distribution  of  a  cannot  be  determined.  However,  the  variance  of  a 
can  be  estimated  by  making  some  simplifying  assumptions. 

Assume  the  X  ,'s  are  normally  and  identically  (but  not  inde- 
J1  2 

pendently)  distributed.  Then  cp..  (0)  is  approximately  X  distributed  with  N 

11  c 

degrees  of  freedom,  where  N  can  be  estimated  from  the  spectral  shape  by 

determining  the  number  of  frequencies  which  contribute  significantly  to  the 

total  power.  Because  of  the  peaked  noise  spectrum,  significant  power  levels 

are  confined  to  about  one-tenth  the  bandwidth,  so  N  ~  -777-  N.  Usually  about 

e  1  u 

2000  points  are  used  to  compute  cp..  (0)  so  that  N  ~  200. 

12.  6 
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Assume  cp^(O)  andcp„(0)  are  independent.  Since  the  reference 

in  this  case  was  seismometer  21,  this  assumption  would  be  approximately 

$ 

true  for  the  seismometers  on  the  outer  edges  of  a  subarray.  Then 

2  V10’ 

a  =  - — 

(0) 

is  approximately  F-distributed  with  200  degrees  of  freedom  for  both  numerator 

1/2 

and  denominator.  Thus ,  a  is  distributed  as  200’ 

The  variance  of  a  is 

Var  (a)  =  E(a2)  -  E2(a) 


■  *H  -  *■[']'" 

■  A  -  *,K‘] 

e 


-  E 


large, 


It  is  shown  in  Appendix  C  that  E 


P/2]~ 


1  when  Ng  becomes 


.  .  Var  (a)  =  N  .  2  1  =  N  -  2 

e  e 


and  the  standard  deviation  of  a  is  approximately  ±10  percent.  Consequently, 
the  assumption  that  a  is  constant  does  not  appear  to  be  reasonable,  and  the 
variations  in  the  equalization  coefficients  cannot  be  attributed  to  seismometer- 
gain  fluctuations  only. 
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Figure  IV- 1  shows  a  plot  of  the  equalization  coefficients  p. 

as  a  function  of  time  obtained  for  several  seismometers.  Note  that  for 

two  noise  samples  on  the  same  day,  considerable  difference  in  o  is  observed, 

i 

implying  that  a  indeed  varies  considerably  from  noise  sample  to  noise  sample. 
Variations  in  p  also  occurred  for  the  inner  seismometers  (where  the  assumption 
that  the  two  time  series  are  independent  does  not  hold),  which  indicates  the 
var.ance  of  a  is  still  significant.  Figure  IV -  1  does  indicate,  however,  that 
the  seismometer  gains  at  low  frequency  are  different  from  channel  to  channel 
because  the  curves  are  separated  by  amounts  larger  than  the  variations 
observed. 

In  summary,  it  is  impossible  to  interpret  variations  in  the 
p/s  solely  in  terms  of  gain  fluctuations.  Instead,  the  variations  are  due  to 
a  combination  of  g;in  fluctuations  and  variations  in  the  estimate  of  a. 


Figure  IV- 1.  Variations  in  Noise  Equalization  Coefficients  for  Seis¬ 
mometers  31,  82,  44,  and  75,  Subarray  B3 
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SECTION  V 

STATISTICAL  PHASE  FLUCTUATIONS 

Considerable  work  was  done  in  the  past  on  the  problem  of 
adding  statistical  gain  fluctuation  to  the  correlation  statistics.  The  purpose 
was  to  prevent  misdesigning  multichannel-filter  systems  by  using  gain  dif¬ 
ferences  to  separate  measured  noise  from  an  equalized  signal  model.  There 
has  been  some  interest  in  the  phase  inequalities  among  sensors  at  the  LASA. 
This  secticr  generalizes  the  method  of  gain  fluctuation  to  incorporate  statisti¬ 
cal  phase  fluctuations  into  the  correlation  matrix. 

Gain  fluctuation  was  incorporated  by  adding  at  each  sensor  a 
1  -  point  stochastic  filter  having  an  expected  value  of  1  and  any  assigned  ari- 
ance.  The  variance  determined  the  amount  of  gain  fluctuation  added.  This 
idea  can  be  generalized  to  gain  and  phase  fluctuation  by  aiding  a  (2L+  l)-point 
stochastic  filter  at  each  seismometer.  These  filters  have  the  following 
properties: 

«•  Expected  value  of  6  (t) 

•  Diagonal  covariance  matrix  (independent 
filter  weights) 

•  No  correlation  between  channels 

To  determine  the  crosscorrelation  between  two  channels 
which  have  stochastic  filters  with  the  proDerties  just  defined,  the  outputs 
of  the  two  channels  are 

f1'(t)=f1  (t)  *  hl  (t) 
f2' (t)  =f2  (t)  *  h2  (t) 
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where 


f  '(t)  and  f  '  (t)  are  outputs  of  channels  1  and  2 

1  2 

h  (t)  and  h  (t)  are  stochastic  filters 

1  Z 

f  (t)  and  f  (t)  are  inputs  to  channels  1  and  2 

1  2 

*  stands  for  convolution 


The  crosscorrelation  is 


$  '  (m)  =  f'  (t)  f/(t  +  m) 

12  L  c 


+L 


f  ^  (t  -  a  A  t)  hj  (o  At) 

a  =  -  L 


2  f  (t  +  m-4At)h2(CAt) 


j£  =  -L 


Writing  out  the  correlation  estimate  as  a  summation,  re¬ 
arranging  terms,  and  taking  the  expectation  inside  the  summation  gives 


$  '(m)  =  ElSf1<nAt-oAt)f2(nAt‘iAt  +  n,|hl(3At)h2^t)  (5_1) 

12  n  l  a 

Since  h^t)  and  h2(t)  are  assumed  independent, 

h  (aAt) h  (£At)  =  h  (crAtj  hzit&t)  =  6  (aAt)  a  (£At) 

1  2  i  - 

and  the  triple  summation  (Equation  5-1)  reduces  to 


*  'M  =  S  fl(nAt)  f2(nAt  +  m)  =  $12(m) 

n 
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Thus,  the  crosscorrelations  are  unchanged  by  the  addition  of  inotpendent 
stochastic  filters  at  each  seismometer. 

The  autocorrelation  is 

$n'{m)  =  f (t)  f  Ct  +  rn) 
and  similarly  can  be  reduced  to  the  form 
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which,  if  end  effects  are  ignored,  reduces  to 


Acl  lu,m) 

L  a  J 


Thus,  all  lags  of  the  original  autocorrelation  are  scaled  by 
the  multiplicative  constant  , 


a  -  -  L 


which  is  equivalent  to  1  plus  the  sum  of  the  variances  assigned  to  the  random 
sample  weights.  This  result  is  the  exact  extension  of  that  obtained  for  the 
1 -point  stochastic  filter,  where  all  lags  of  the  autocorrelation  were  scaled 
by  1  plus  the  variance  of  the  random  sample  weight. 
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appendix  a 

DESCRIPTION  OF  PROGRAM  TO  COMPUTE 
MICROSTATIC  CORRECTIONS  ON  SIGNAL  TRACES 

This  program  is  written  so  that  the  signal  traces  can  be 
aligned  to  a  fraction  of  a  sample,  which  is  necessary  in  order  to  use  the 
difference  between  the  reference  and  signal  traces  as  a  method  of  evalu¬ 
ating  the  effect  of  equalization.  The  program  is  divided  into  three  steps. 

(1)  The  traces  are  statically  aligned  to  the  nearest 
integer  sample 

(2)  The  best  estimate  of  the  microstatic  shift  a 
is  determined 

(3)  The  interpolation  filter  which  applies  the 
shift  a  is  designed  and  applied  to  the  traces 

The  maximum  lag  of  the  signal-reference  crosscorrelation 
is  used  to  align  the  signals  to  the  nearest  integer  sample. 

To  estimate  a,  treat  the  crosscorrelation  function  between 
the  signal  and  reference  traces  cpgT  as  an  ordinary  time  trace.  Then  find 
the  Wiener  alter  which  estimates  the  crosscorrelation  at  T  =  (n+  1/2)  At. 
That  is,  solve 

!>][*]=  W 

where 

cp  is  autocorrelation  of  the  crosscorrelation 
f  is  interpolation  filter 

c  is  inverse  transform  of  the  "spectrum"  of 
the  crosscorrelation  at  T  =  (n+  1/2)  At, 
n  =  0,  ±  1  •  •  • 
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It  can  be  shown  that  cp  can  be  obtained  by  convolving  the  autocorrelation  of 
the  signal  and  target  traces;  i.  e.  , 


cp  =  [X(t)  cp  T(t)]  cp  [X(t)  cp  T  (t)] 
=  [x(t)cpX(t)]  *  [t  (t)  cp  T  (t)] 


where 

X  (t)  is  signal  trace 
T  (t)  is  target  trace 
cp  is  correlation 
*  is  convolution 

The  interpolation  filter  then  is  applied  to  tPgT>  and  the  value 
of  the  maximum  crosscorrelation  lag  tp(±-y-)  is  compared  with  cp(0).  Assume 
cp^-^-^is  larger  than  cp  y~ )•  Then,  if  the  inequality 

|cp(0)  -  cp(y  At)  | 

_ -  <  10 

max[cp(0), 

is  satisfied,  a  is  chosen  as  either  0  or  At/2,  whichever  has  the  largest 
crosscorrelation  value.  Otherwise,  the  process  is  repeated  until  the  in¬ 
equality  is  satisfied  or  until  15  iterations  are  completed  (usually  six  or 
seven  iterations  are  sufficient). 

Having  chosen  a,  we  proceed  to  the  third  step,  which  is  to 
compute  the  interpolation  filter  for  the  signal  trace.  This  is  achieved  by 

solving 

[•JM  ■  M 


A-2 


•clones  •srvlcsa  division 


Mi 


where 


cp  is  autocorrelation  of  X(t) 

xx 

T  is  interpolation  filter 

K  is  inverse  transform  of  the  spectrum 

X(t)  at  t  =  nAt  +  a,  n=0,  ±1... 


This  is  exactly  the  same  equation  as  in  the  second  step,  except  here  we  are 
dealing  with  the  signal  trace  instead  of  the  crosscorrelation  between  the 
signal  and  target  traces.  The  filter  T  then  is  applied  to  X(t)  to  microstatically 
align  it  with  the  reference  trace. 
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APPENDIX  B 

PROBABILITY  DISTRIBUTION 
OF  THE  RATIO  OF  INDEPENDENT, 
IDENTICALLY  DISTRIBUTED,  UNIFORM, 
RANDOM  VARIABLES 


Suppose  X  and  Y  are  independent,  identically  distributed, 
random  variables  with  probability  density  function  (PDF) 


f  (x) 


a  <  x  <  b 
a  <  y  <  b 


=  V 


elsewhere 


We  want  to  find  the  PDF  of  Z  =  X/Y.  Let  W  =  X  and  Z  =  X/Y .  Solving  for 
X  and  Y  in  terms  of  W  and  Z  gives  X  =  W  and  Y  =  W/Z.  Then, 


'  ax 

BX  “ 

1 

0 

aw 

BZ 

BY 

BY 

1 

w 

_aw 

BZ  _ 

_  Z 

z  . 

Therefore, 


f  (w,  z) 


1 


w 


2 


z 


The  mapping  from  the  XY  to  WZ  plane  is  as  follows. 
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To  find  f(z),  integrate  out  w  in  f(w,  z)  as  follows: 


f  (z) 


•  • 


f  (z) 


1 

2 

2 

(b  -  a) 

z 

1 

,2 

2 

(b  -  a) 

z 

1 

2 

2 

(b  -  a) 

z 

1 

2  (b  -  a)2  z2 


w  dw 


w  dw 


a 

b 


<  z  <  1 


1  <  z  < 


b 

a 


a_ 

b 


<  z  <  1 


_1_ 

2 


(b  -  a) 


(b 


.2  2 

a)  z 


0 

0 

0 

0 

0 

0 

[] 


Similarly, 


f  (z)  = 


■+[-, 


(b  -  a)2  z2  (b  -  a)2 


] 


Therefore,  the  PDF  of  z  is 

.  2 


f  (z) 


_1 

2 


1 


2  ,  ,2  2  , 

(b  -  a)  (b  -  a)  z  J 


(b  -  a)2  z2  (b  -  a) 


] 


_1_ 

2 


=  0 


Check  to  determine  that  PDF  integrates  to  1, 


1  <  z  <  — 
a 


—  <  z  <  1 


1  <  z  <  — 
a 


elsewhere 


u 

0 

0 


0 


u 


f £(a| dz 


dz 


a/b  (b  -  a) 


•I. 


2  2 

a/b  (b  -  a)  z 


dz 


,b/a  ,  2  /*b/a  2 

+  I  - rr  dz  -  J  — “ — “  dz 

•'z  =  1 


/•Old 

X  =  i 


2  2 

(b  -  a)  z 


(b  -  a) 


■+[■ 


a  b  a 


b^a  b-a  b-a  b 


h] 


1  fb  -  a  b  -  a  ~l  _  _1_ 

2  |_b  -  a  b-aj  2 


2  =  1 
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APPENDIX  C 

EXPECTED  VALUE  OF  THE  SQUARE  ROOT  OF  F 


We  want  to  find  the  expected  value  of  F  ,  where  F  is 
F-distributed  with  n  degrees  of  freedom  for  both  numerator  and  denomi¬ 
nator.  For  simplicity,  assume  that  n  is  even  (as  it  was  for  all  compu¬ 
tations).  The  PDF  of  F  is 


n- 2 


The  expression  in  braces  is  the  integral  of  an  F-distribution 
with  n  +  1,  n  -  1  degrees  of  freedom  over  the  whole  range  of  the  random 
variable  and  must  be  unity. 


Using  Stirling's  approximation  for  n!, 
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